Multisample analysis#

Setup#

import missionbio.mosaic as ms
from IPython.display import display, HTML

# Load a multi-sample file
# By default, a multi-sample h5 file loads as a SampleGroup object

group = ms.load_example_dataset("Multisample PBMC")

# To analyze one sample at a time access, them using the samples attribute
sample = group.samples[1]

# These are the samples in the h5 file
[s.name for s in group.samples]

['Sample 1', 'Sample 2']

Applying functions#

# To apply a function across all samples use the `apply` method on the SampleGroup
# It returns a list of returned objects for each sample

def filt(sample):
    filt_vars = sample.dna.filter_variants()
    return filt_vars

filtered_variants = group.apply(filt)
filtered_variants

[array(['chr2:198266943:C/T', 'chr2:198267770:G/GAA', 'chr4:55599436:T/C',
        'chr5:170837457:A/G', 'chr7:140449071:C/G', 'chr7:148504716:AG/A',
        'chr7:148504854:A/AGACTT', 'chr7:148508833:A/G',
        'chr7:148543525:A/G', 'chr7:148543583:G/C', 'chr7:148543693:TA/T',
        'chr11:32417945:T/C', 'chr11:119148573:G/T', 'chr13:28602292:T/C',
        'chr17:7578115:T/C', 'chr17:7579801:G/C', 'chr17:29559932:C/A',
        'chr20:31024028:G/A', 'chrX:39933339:A/G', 'chrX:44833841:C/A',
        'chrX:133547814:T/C'], dtype='<U23'),
 array(['chr2:198266943:C/T', 'chr2:198267770:G/GAA', 'chr4:55599436:T/C',
        'chr5:170837457:A/G', 'chr7:140449071:C/G', 'chr7:148504716:AG/A',
        'chr7:148504854:A/AGACTT', 'chr7:148508833:A/G',
        'chr7:148543525:A/G', 'chr7:148543583:G/C', 'chr7:148543693:TA/T',
        'chr11:32417945:T/C', 'chr11:119148573:G/T', 'chr13:28602292:T/C',
        'chr17:7578115:T/C', 'chr17:7579801:G/C', 'chr17:29559932:C/A',
        'chr20:31024028:G/A', 'chrX:39933339:A/G', 'chrX:44833841:C/A',
        'chrX:133547814:T/C'], dtype='<U23')]

# Subset the same variants in all dna assays
# It is important to maintain the same variants across all dna assays

og_num_vars = [s.dna.shape[1] for s in group.samples]

var_union = list(set().union(*filtered_variants))
for sample in group:
    sample.dna = sample.dna[:, var_union]  # Subsets all samples with the same variants

new_num_vars = [s.dna.shape[1] for s in group.samples]

print(og_num_vars, new_num_vars)  # Thee old and new number of variants for each sample in the group

[21, 21] [21, 21]

# The functions applied on each sample can be more complex - like this assignment and relabeling method
# Note the original labels can be uncoordinated across samples in the group
# The labels are changed to ensure that each label is for the same clone

variants_of_interest = ['chr7:148508833:A/G', 'chr17:29559932:C/A', 'chr4:55599436:T/C']
def cluster(sample):
    clone_table = sample.dna.group_by_genotype(variants_of_interest, max_ado_score=0.8)

    # Rename labels so that each sample has the same labels
    # Here the signature of each variant is used to rename the labels
    df = sample.dna.signature("NGT").loc[:, variants_of_interest]
    names = df.apply(lambda vs: "-".join([str(int(v)) for v in vs]), axis=1)
    label_map = {i: n for i, n in names.items()}
    
    # Don't rename the outlier categories:
    for lab in ["missing", "small", "ADO"]:
        del label_map[lab]

    sample.dna.rename_labels(label_map)
    clone_table = clone_table.rename(index=label_map)

    return clone_table  # Return the clone table

tables = group.apply(cluster)

for t in tables:  # The clone tables for each sample
    display(HTML(t.to_html()))
clone 1 3 Missing GT clones (27) Small subclones (10) ADO clones (4)
chr7:148508833:A/G Het (54.97%) WT (0.81%) Missing in 20.64% of cells Het (43.28%) Het (52.45%)
chr17:29559932:C/A Het (55.33%) Het (58.88%) Missing in 19.25% of cells Het (57.55%) Het (55.48%)
chr4:55599436:T/C Hom (99.88%) Het (56.63%) Missing in 19.71% of cells Het (75.64%) Hom (99.86%)
Total Cell Number 769 (39.58%) 42 (2.16%) 929 (47.81%) 43 (2.21%) 160 (8.23%)
Sample 1 Cell Number 769 (39.58%) 42 (2.16%) 929 (47.81%) 43 (2.21%) 160 (8.23%)
Parents NaN NaN NaN NaN NaN
Sisters NaN NaN NaN NaN NaN
ADO score 0 0 NaN NaN NaN
clone 1 3 5 Missing GT clones (24) Small subclones (7) ADO clones (4)
chr7:148508833:A/G WT (0.79%) Het (53.93%) Het (39.14%) Missing in 8.35% of cells WT (24.64%) WT (0.63%)
chr17:29559932:C/A Het (54.45%) Het (55.08%) Het (53.12%) Missing in 19.69% of cells (60.04%) Het (57.45%)
chr4:55599436:T/C Het (55.97%) Hom (99.89%) Het (71.68%) Missing in 35.23% of cells Hom (74.7%) Het (55.76%)
Total Cell Number 646 (35.73%) 38 (2.1%) 29 (1.6%) 945 (52.27%) 24 (1.33%) 126 (6.97%)
Sample 2 Cell Number 646 (35.73%) 38 (2.1%) 29 (1.6%) 945 (52.27%) 24 (1.33%) 126 (6.97%)
Parents NaN NaN NaN NaN NaN NaN
Sisters NaN NaN NaN NaN NaN NaN
ADO score 0 0 0 NaN NaN NaN

Drawing figures#

# Displaying the same plot across all the samples

for sample in group:
    sample.dna.heatmap("NGT_FILTERED", features=variants_of_interest).show()

# Normalize Protein and look for patterns

for sample in group:
    sample.protein.normalize_reads()
    sample.heatmap(("dna", "protein")).show()

variants_of_interest = ['chr17:29559932:C/A', 'chr4:55599436:T/C', 'chr7:148508833:A/G']
proteins_of_interest = ['CD19', "CD34", "CD30"]
clones_of_interest = ["0-1-1", "1-1-2"]

for sample in group:
    s = sample[sample.dna.barcodes(clones_of_interest)]
    s.dna = s.dna[:, variants_of_interest]
    s.protein = s.protein[:, proteins_of_interest]
    s.clone_vs_analyte("protein")
../_images/ca58f51dc1f2cd6b249efdffc9e24afd938e785fea673845c2efc7c127c62ee7.png ../_images/2dce9837bedb5f538242dccec61c9f0aec12cd464cacac84b8db87b686e9121f.png

Multisample plots#

# Draw a fishplot for the dna labels
# From the proportions in the heatmaps - The two clones of interest are 0-1-1 and 2-20

group.fishplot(
    "dna",
    sample_order=["Sample 1", "Sample 2"],
    labels=["0-1-1", "1-1-2"],
    parents=[None, None]
)

# Draw a barplot for the dna labels

group.barplot(
    "dna",
    sample_order=["Sample 1", "Sample 2"],
    label_order=["0-1-1", "1-1-2"],
    percentage=True
)